Matrix exponential and Krylov subspaces for fast time domain computations: recent advances

We show how finite difference (or finite element) time domain computations can be accelerated by employing recent advances in the matrix exponential time integration and Krylov subspace techniques

A block Krylov subspace time-exact solution method for linear ODE systems

We propose a time-exact Krylov-subspace-based method for solving linear ODE (ordinary differential equation) systems of the form $y'=-Ay + g(t)$ and $y''=-Ay + g(t)$, where $y(t)$ is the unknown function. The method consists of two stages. The first stage is an accurate piecewise polynomial approximation of the source term $g(t)$, constructed with the help of the truncated SVD (singular value decomposition). The second stage is a special residual-based block Krylov subspace method. The accuracy of the method is only restricted by the accuracy of the piecewise polynomial approximation and by the error of the block Krylov process. Since both errors can, in principle, be made arbitrarily small, this yields, at some costs, a time-exact method. Numerical experiments are presented to demonstrate efficiency of the new method, as compared to an exponential time integrator with Krylov subspace matrix function evaluations

Residual, restarting and Richardson iteration for the matrix exponential

A well-known problem in computing some matrix functions iteratively is a lack of a clear, commonly accepted residual notion. An important matrix function for which this is the case is the matrix exponential. Assume, the matrix exponential of a given matrix times a given vector has to be computed. We interpret the sought after vector as a value of a vector function satisfying the linear system of ordinary differential equations (ODE), whose coefficients form the given matrix. The residual is then defined with respect to the initial-value problem for this ODE system. The residual introduced in this way can be seen as a backward error. We show how the residual can efficiently be computed within several iterative methods for the matrix exponential. This completely resolves the question of reliable stopping criteria for these methods. Furthermore, we show that the residual concept can be used to construct new residual-based iterative methods. In particular, a variant of the Richardson method for the new residual appears to provide an efficient way to restart Krylov subspace methods for evaluating the matrix exponential.\u

Experiments with MRAI time stepping schemes on a distributed memory parallel environment

Implicit time stepping is often difficult to parallelize. The recently proposed Minimal Residual Approximate Implicit (MRAI) schemes are specially designed as a cheaper and parallelizable alternative for implicit time stepping. A several GMRES iterations are performed to solve approximately the implicit scheme of interest, and the step size is adjusted to guarantee stability. A natural way to apply the approach is to modify a given implicit scheme in which one is interested. Here, we present numerical results for two parallel implementations of MRAI schemes. One is based on the simple Euler Backward scheme, and the other is the MRAI-modified multistep ODE solver LSODE. On the Cray T3E and IBM SP2 platforms, the MRAI codes exhibit parallelism of explicit schemes. The model problem under consideration is the 3D spatially discretized heat equation. The speed-up results for the Cray T3E and IBM SP2 are reported

Improving approximate matrix factorizations for implicit time integration in air pollution modelling

For a long time operator splitting was the only computationally feasible way of implicit time integration in large scale Air Pollution Models. A recently proposed attractive alternative is Rosenbrock schemes combined with Approximate Matrix Factorization (AMF). With AMF, linear systems arising in implicit time stepping are solved approximately in such a way that the overall computational costs per time step are not higher than those of splitting methods. We propose and discuss two new variants of AMF. The first one is aimed at yet a further reduction of costs as compared with conventional AMF. The second variant of AMF provides in certain circumstances a better approximation to the inverse of the linear system matrix than standard AMF and requires the same computational work

On Skew-symmetric Preconditioning for Strongly Non-symmetric Linear Systems

To solve iteratively linear system $Au=b$ with large sparse strongly non-symmetric matrix $A$ we propose preconditioning $\hat A \hat u = \hat b$, $\hat A=(I+\tau L_1)^{-1} A (I+\tau U_1)^{-1},\; \tau>0$ where respectively lower and upper triangular matrices $L_1$ and $U_1$ are so that $L_1+U_1=1/2(A-A^*)$. Such preconditioning technique may be treated as a variant of ILU-factorization, and we call it MSSILU --- Modified Skew-Symmetric ILU. \ud We investigate and optimize (with respect to $\tau$) convergence of preconditioned Richardson method (RM) of the following special form: ${\hat x}^{m+1}=(I-\tau \hat A){\hat x}^m+\tau {\hat b},\; m\geq 0$, where $\tau$ is the same as in $\hat A$. For this method we give an estimate for rate of convergence in relevant Euclidean norm for the case of positivereal matrix $A$. \ud Numerical experiments have included solving linear systems arising from 5-point FD approximation of convection--diffusion equation with dominated convection by MSSILU+RM, MSSILU+GMRES(2) and MSSILU+GMRES(10).\u

Parallel processing and non-uniform grids in global air quality modeling

A large-scale global air quality model, running efficiently on a single vector processor, is enhanced to make more realistic and more long-term simulations feasible. Two strategies are combined: non-uniform grids and parallel processing. The communication through the hierarchy of non-uniform grids interferes with the inter-processor communication. We discuss load balance in the decomposition of the domain, I/O, and inter-processor communication. A model shows that the communication overhead for both techniques is very low, whence non-uniform grids allow for large speed-ups and high speed-up can be expected from parallelization. The implementation is in progress, and results of experiments will be reported elsewhere

A new approximate matrix factorization for implicit time integration in air pollution modeling

Implicit time stepping typically requires solution of one or several linear systems with a matrix I−τJ per time step where J is the Jacobian matrix. If solution of these systems is expensive, replacing I−τJ with its approximate matrix factorization (AMF) (I−τR)(I−τV), R+V=J, often leads to a good compromise between stability and accuracy of the time integration on the one hand and its efficiency on the other hand. For example, in air pollution modeling, AMF has been successfully used in the framework of Rosenbrock schemes. The standard AMF gives an approximation to I−τJ with the error τ2RV, which can be significant in norm. In this paper we propose a new AMF. In assumption that −V is an M-matrix, the error of the new AMF can be shown to have an upper bound τ||R||, while still being asymptotically $O(\tau^2)$. This new AMF, called AMF+, is equal in costs to standard AMF and, as both analysis and numerical experiments reveal, provides a better accuracy. We also report on our experience with another, cheaper AMF and with AMF-preconditioned GMRES

A class of nonsymmetric preconditioners for saddle point problems

For iterative solution of saddle point problems, a nonsymmetric preconditioning is studied which, with respect to the upper-left block of the system matrix, can be seen as a variant of SSOR. An idealized situation where the SSOR is taken with respect to the skew-symmetric part plus the diagonal part of the upper-left block is analyzed in detail. Since action of the preconditioner involves solution of a Schur complement system, an inexact form of the preconditioner can be of interest. This results in an inner-outer iterative process. Numerical experiments with solution of linearized Navier-Stokes equations demonstrate efficiency of the new preconditioner, especially when the left-upper block is far from symmetric